Data Citations: Land use/land cover data: http://geoportal.hawaii.gov/datasets/land-use-land-cover-lulc Watershed data: http://geoportal.hawaii.gov/datasets/watersheds All other hydrological data:http://geoportal.hawaii.gov/datasets
## Attach packages:
library(tidyverse)
library(janitor)
library(lubridate)
library(here)
library(paletteer)
library(sf)
library(tmap)
library(mapview)
library(tmaptools)
library(rmapshaper)
Read in original shapefiles from Assignment 2 - Part 3 and ensure all are in same crs (coordinate system and projection).
Watersheds<- read_sf(dsn = here::here("HI Data", "hw watersheds"),
layer = "Watersheds") %>%
st_transform(crs = 4326) %>%
select("wuname", "area_sqmi")
land <- read_sf(dsn = here::here("HI Data","Land_Use_Land_Cover_LULC"),
layer = "Land_Use_Land_Cover_LULC") %>%
st_transform(crs = 4326) %>%
filter(landcover != 0)
Water <- land %>% #Making a dataframe that only includes water oriented landcovers.
filter(landcover == "Streams and Canals" | landcover == "Forested Wetland" | landcover == "Resevoirs" | landcover == "Nonforested Wetland" | landcover == "Lakes" | landcover == "Bays and Estuaries") %>%
rename(Landcover = landcover)
#Exploring new water dataframe:
#mapview(Water)
#Interestingly, not much is shown on each island for the different landcovers under "Water."
Industry <- land %>% #Making a dataframe that only includes man-made/industrial/pollution oriented landcovers.
filter(landcover == "Residential" | landcover == "Other Urban or Built-up Land" | landcover == "Industrial" | landcover == "Industrial and Commercial Complexes" | landcover == "Commercial and Service" | landcover == "Transportation, Communications and Utilities" | landcover == "Other Agricultural Land" | landcover == "Strip Mines, Quarries, and Gravel Pits" | landcover == "Mixed Urban or Built-up Land" | landcover == "Confined Feeding Operations")
#Exploring new industry datafram:
#mapview(Industry)
#Similar to "Water", not much is shown. While industral landcovers are where I would think (near city centers, etc.), there is clearly a lot missing.
Read in additional shapefiles to explore hydrology of Hawaii, to contrast with water oriented landcover dataframe, and ensure all are in same crs (coordinate system and projection).
dlnr_aquifers<- read_sf(dsn = here::here("HI Data","DLNR_Aquifers_Poly"), #More resource-oriented in nature aquifer boundaries.
layer = "DLNR_Aquifers_Poly") %>%
st_transform(crs = 4326)
doh_aquifers <- read_sf(dsn = here::here("HI Data","DOH_Aquifers_Polygons"), #More administrative in boundaries. Interesting that they have these two different datasets.
layer = "DOH_Aquifers_Polygons") %>%
st_transform(crs = 4326)
nhd_flowlines <- read_sf(dsn = here::here("HI Data","NHD_Flowlines"), #All HUC 8 flowlines.
layer = "NHD_Flowlines") %>%
st_transform(crs = 4326)
ms_simplify(nhd_flowlines) #Simplifying flowlines because working with this dataset is too difficult/time consuming as is.
## Simple feature collection with 16954 features and 16 fields (with 2268 geometries empty)
## geometry type: LINESTRING
## dimension: XY
## bbox: xmin: -160.2466 ymin: 18.91069 xmax: -154.8069 ymax: 22.23272
## epsg (SRID): 4326
## proj4string: +proj=longlat +datum=WGS84 +no_defs
## # A tibble: 16,954 x 17
## objectid permanent_ fdate resolution gnis_id gnis_name lengthkm reachcode
## <dbl> <chr> <chr> <dbl> <chr> <chr> <dbl> <chr>
## 1 1 31079000 2012… 2 <NA> <NA> 5.33 20070000…
## 2 2 31079692 2012… 2 <NA> <NA> 5.39 20070000…
## 3 3 31079996 2013… 2 <NA> <NA> 4.76 20070000…
## 4 4 31086060 2013… 2 <NA> <NA> 7.56 20070000…
## 5 5 31082174 2016… 2 003618… Lāwaʻi S… 4.26 20070000…
## 6 6 31071620 2012… 2 <NA> <NA> 5.59 20070000…
## 7 7 28026375 2012… 2 <NA> <NA> 6.44 20030000…
## 8 8 32453113 2012… 2 <NA> <NA> 14.8 20020000…
## 9 9 32447113 2013… 2 <NA> <NA> 5.59 20020000…
## 10 10 32447079 2013… 2 <NA> <NA> 7.29 20020000…
## # … with 16,944 more rows, and 9 more variables: flowdir <dbl>,
## # wbarea_per <chr>, ftype <dbl>, fcode <dbl>, mainpath <dbl>,
## # innetwork <dbl>, visibility <dbl>, st_lengths <dbl>, geometry <LINESTRING
## # [°]>
#mapview(nhd_flowlines)
Streams <- read_sf(dsn = here::here("HI Data","Streams"),
layer = "Streams") %>%
st_transform(crs = 4326) %>%
mutate(Stream = str_to_lower(type)) #Type: Perennial (continuous flow) and Non-Perennial (seasonal flows).
Rain <- read_sf(dsn = here::here("HI Data","StateIsohyetsSHP_mm"), #Rain contour lines for all islands, might be interesting to contrast against streams (perennial vs. non).
layer = "isohyet_mm_01") %>%
st_transform(crs = 4326)
Create clips of a couple dataframes for more convinient ggplots (geom_sf) later.
#Clips to explore in mapping below.
Flowlines <- st_intersection(Watersheds, nhd_flowlines) %>% #This was a very large file (long time to run), however it works more efficiently with simplification above.
select("wuname", "gnis_name", "st_lengths")
Waterways <- st_intersection(Watersheds, Water) %>%
select("wuname", "Landcover")
kauai_watercover <- ggplot(data = Watersheds) +
geom_sf(color = "grey3",
size = 0.2) +
geom_sf(data = Water,
aes(fill = Landcover),
show.legend = TRUE,
color = "NA",
alpha = 0.8) +
ggtitle("Water Landcovers on Kauai, HI") +
coord_sf(xlim = c(-159.2, -159.9), ylim = c(21.8, 22.3), expand = FALSE) + #Cropping Kauai as example island to see what is going on.
theme_classic()
kauai_watercover
kauai_aquifers <- ggplot() +
geom_sf(data = doh_aquifers,
aes(colour = "DOH (Administratively Oriented)"),
alpha = 0.5,
show.legend = "line") +
geom_sf(data = dlnr_aquifers,
aes(colour = "DLNR (Hydrologically Oriented)"),
alpha = 0.5,
show.legend = "line") +
coord_sf(xlim = c(-159.2, -159.9), ylim = c(21.8, 22.3), expand = FALSE) + #Cropping Kauai as example island to see what is going on.
ggtitle("Aquifer Boundary Deliniations on Kauai, HI\n ") +
scale_colour_manual(values = c("DLNR (Hydrologically Oriented)" = "palegreen4", "DOH (Administratively Oriented)" = "tomato2")) +
guides(color=guide_legend("Tyep of Aquifer Boundary")) +
theme_classic()
kauai_aquifers
kauai_streams <- ggplot(data = Flowlines) +
geom_sf(aes(colour = wuname),
alpha = 0.5,
show.legend = FALSE) +
ggtitle("Flowlines by Watershed on Kauai, HI")+
theme_classic() +
coord_sf(xlim = c(-159.2, -159.9), ylim = c(21.8, 22.3),expand = FALSE) #Cropping Kauai as example island to see what is going on.
kauai_streams
#Included watersheds to determine if larger watershed have more perennial streams (only based on visual cues). This is not possible however because the data seems questionable. On many of the islands, very small watersheds are shown to have much more area than larger. The square mileage was computed with GIS, and perhaps there was an issue there.
tmap_mode("view") #Set so the map will be interactive.
all_islands_tmap <- tm_basemap("Esri.WorldImagery") + #Add a basemap explored previously in mapview.
#tm_shape(Watersheds) +
#tm_polygons("area_sqmi",
# style = "pretty",
# palette = "YlOrRd",
# id = "wuname",
# popup.var = TRUE,
# colorNA = NULL,
# title = "Watershed Area (sq. miles)",
# lwd = 1,
# alpha = 0.5) +
tm_shape(Streams) +
tm_lines(col = "Stream",
title.col = "Stream Type",
style = "pretty",
palette = "Greens",
id = "Stream",
popup.var = FALSE) +
tm_shape(Rain) +
tm_lines(col = "CONTOUR",
title.col="Rain Contour (mm)",
style = "pretty",
palette = "Blues",
id = "Contour",
popup.var = TRUE,
colorNA = NULL) +
tm_layout("Hawaiian Island Hyrdological Dynamics (WGS84)")
all_islands_tmap